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IN-PLANE  DEFORMATION  OF 
NON-COAXIAL  PLASTIC  SOIL 

Shunsuke  Takagi 


INTRODUCTION 

The  work  on  mathematical  soil  plasticity  started  with  Drucker  and  Prager.’  Shield,”  following 
their  theory,  showed  some  typical  solutions  of  the  in-plane  deformation.  However,  the  plastic  de- 
formation defined  by  these  authors  is  always  accompanied  by  dilatancy,  and  deviates  from  the  slip 
direction  by  a definite  angle.  The  solution  of  the  deformation  derived  by  them  does  not  satisfy  the 
boundary  slip-line  condition  that  imposes  restrictions  on  both  stress  and  flow. 

Haythornthwaite,”  Dejosselin  de  Jong,*  and  Mandl”  expressed  suspicion  and  suggested  ways 
of  modifying  Drucker  and  Prager’s  implicit  assumption  that  the  stress  and  strain  rate  tensors  are 
coaxial.  If  the  coaxiality  is  no  longer  true,  the  strain-rate  tensor  is  not  normal  to  the  yield  criterion 
surface  in  the  six-dimensional  stress  space,  and  the  theory  of  plastic  potential  (Drucker;*  Drucker,* 
p.  273)  cannot  bt  used  to  derive  the  constitutive  equations.  Then,  the  basic  premise  of  Drucker  and 
Prager’s’  developn^^ent  does  not  hold. 

On  the  other  hand,  Geniev*  formulated  a theory  of  in-plane  deformation  assuming  that  the  non- 
coaxial angle  is  equa'  to  one  half  the  frictional  angle.  The  deformation  he  formulated  causes  no 
volume  change.  We“  extended  the  non<oaxial  angle  to  an  arbitrary  value.  The  deformation  we 
formulated  can  accommodate  both  dilatancy  and  compression.  However,  introduction  of  the  un- 
known non-coaxial  angle  increases  the  number  of  unknowns  by  one;  no  method  of  evaluating  the 
non-coaxial  angle  has  ever  been  discovered,  and  no  bounr*  iry-value  problem  has  ever  been  solved. 
This  difficulty  is  now  overcome.  We  discovered  that  the  non-coaxial  angle  must  be  so  determined 
that  the  presupposed  flow  can  be  accommodated  in  the  given  domain  satisfying  the  given  boundary 
conditions. 

In  the  following,  first  we  introduce  the  geometry  of  the  Mohr  circle,  which  is  developed  into  a 
complete  mathematical  tool  compatible  with  the  analytical  method.  This  is  achieved  by  introducing 
the  sign  of  tangential  components,  called  sign  m,  which  we  believe  will  find  applications  also  in 
other  branches  of  mechanics  than  in  soil  mechanics. 

Second,  we  derive  the  constitutive  equations  by  applying  the  geometry  of  the  Mohr  circle  and 
the  theory  of  the  characteristic  lines.  In  our  theory,  it  is  not  required  that  both  of  the  strain-rate 
characteristic  directions  be  coincident  with  both  of  the  stress  characteristic  directions.  For  a plastic 
deformation  to  occur,  it  is  sufficient  that  one  member  in  a set  of  strain-rate  characteristic  directions 
be  coincident  with  the  sane  sign  m member  of  the  stress  characteristic  directions.  The  coincident 
and  noncoincident  directions  are  called  doublet  and  singlet,  respectively.  In  this  way,  we  formulate 
the  volume  characteristic  equation,  which  shows  that  either  volume  expansion,  shrinkage,  or  no 
volume  change  can  occur  in  accordance  with  the  value  of  the  non-coaxial  angle.  Introduction  of  this 
equation  enables  us  to  complete  the  formulation  of  the  constitutive  equations. 

The  differential  equations  for  solving  a boundary  value  problem  in  this  paper  are  derived  by  use 
of  the  following  property  of  the  strain-rate  characteristic  lines:  the  lengths  of  the  curve  elements  of 
the  strain-rate  characteristic  lines  are  maintained  constant  during  the  plastic  deformation.  The 
expression  of  this  fact,  or  more  in  general,  the  expression  of  the  components  of  the  simplified 


strain-rate  tensor,  by  use  of  an  appropriately  chosen  coordinate  system  gives  the  differential  equa- 
tions for  solving  a deforwjtion  boundary  value  problem.  The  differential  equations  contain  the 
non-coaxial  angle,  to  which  we  can  assign  an  appropriate  value  so  that  the  presupposed  type  of  flow 
can  be  fitted  into  the  given  domain.  As  an  example,  Hill’s"  solution  of  the  stress  distribution  in 
the  ground  sustaining  a rectangular  load  is  given  a deformation  solution.  The  solution  presented 
here  is  the  simplest  one  that  can  be  fitted  into  the  given  domain  of  the  stress  solution. 

We  introduce  the  term  Coulomb  material  to  designate  the  non-coaxial  plastic  material  that  obeys 
the  Coulomb  yield  criterion  in  the  in-plane  deformation.  The  concept  of  Coulomb  material  enables 
us,  we  believe,  to  start  the  rational  mechanics  formulation  of  soil  deformation.  However,  before  we 
can  discover  a satisfactory  model  of  soil  deformation,  we  must  endow  the  Coulomb  material  with 
some  more  rational  mechanics  concepts,  as  discussed  in  the  conclusion  of  this  report. 

In  this  report,  we  have  revised  several  concepts,  improved  the  terminology,  and  corrected  mis- 
takes in  our  previous  papers.”  “ ” 


ANALYSIS  OF  STRESS 
Geometry  of  the  Mohr  circle 

Terzaghi^  introduced  the  concept  of  pole  on  the  Mohr  circle  to  describe  the  orientation  of  a 
stress  tensor  in  the  physical  plane.  In  Figure  la,  A and  B are  the  points  expressing  the  major  and 
minor  principal  stresses,  o,  and  O2,  respectively.  In  the  physical  plane,  there  exists  a pair  of  lines  on 
which  o,  and  Oj  work.  (They  are  shown  by  A,  Q,  and  B,  0]  on  the  right  side  of  Figure  lb,  and 
by  Aj  R2  and  Bj  R2  on  the  left  side,  respectively.)  Let  us  draw  lines  AP  and  BP  parallel  to  the  lines 
on  whose  normal  directions  Oj  and  Oj,  respectively,  work.  The  lines  AP  and  BP  intersect  at  point 
P on  the  Mohr  circle,  which  is  called  pole.  Pole  has  the  following  property  (Terzaghi,**  p.  1 8):  Let 
Q be  a point  on  the  Mohr  circle  representing  stress  (o,  r).  Then,  line  PQ  is  parallel  to  a physical  line 
on  which  the  stress  (a,  t)  in  the  physical  plane  works. 

Locating  the  stress  point  on  either  the  upper  or  the  lower  semicircle  is  essentially  important  to 
the  geometry  of  the  Mohr  circle.  In  the  following,  an  effective  method  that  enables  us  to  decide 
the  location  of  the  stress  point  will  be  presented. 

First,  we  consider  a point  Q on  the  upper  semicircle  in  Figure  la.  We  draw,  as  shown  in  the  right 
side  of  Figure  1b,  0,  A, , Q,  B, , and  A,  B,  parallel  to  PA,  PB,  and  PQ.  Note  that  the  hypotenuse 
A)  B,  is  parallel  to  the  physical  line  under  consideration.  The  tangential  stress  t on  both  sides  of 
A,  B]  makes  the  rotational  direction  as  shown,  as  may  be  verified  by  decomposing  o,  and  02  in  the 
directions  parable'  and  normal  to  A,  Bj . To  designate  the  sign  of  t,  we  suppose  that  a pair  of  r's  on 
both  sides  o^  A,  B,  form  a couple.  We  find  that  the  fictitious  couple  formed  by  the  tangential  stress 
components  makes  the  counterclockwise  rotational  direction. 

Second,  we  consider  point  R on  the  lower  semicircle  of  Figure  la.  Repeating  the  similar  pro- 
cedure, as  shown  on  the  left  of  Figure  1b,  we  find  that  the  fictitious  moment  formed  by  a pair  of 
r’s  on  both  sides  of  A2  B2,  A2  B2  being  parallel  to  PR,  makes  the  clockwise  rotational  direction. 

The  direction  of  r thus  defined  is  convenient  for  locating  a stress  point  on  the  Mohr  circle.  Let 
us  suppose  a fictitious  tangential-stress  moment  formed  by  the  slipping  motion  on  a slip-line  bound- 
ary; then,  we  can  locate  the  stress  point  of  a slip-line  boundary  on  the  Mohr  circle  (see  Fig.  6a  and 
6b).  This  clue  usually  enables  us,  even  prior  to  solving  the  problem,  to  locate  the  stress  points  of  the 
entire  region  on  the  Mohr  circle. 

The  counterclockwise  and  clockwise  rotational  directions  of  the  fictitious  moment  are  called 
siffn  +m  and  -m,  respectively.”  “ The  +/»  and  -m  stress  points  are  on  the  upper  and  lower  semi- 
circles, respectively,  of  the  stress  Mohr  circles. 

The  geometry  of  the  Mohr  circle  thus  amplified  with  the  introduction  of  the  sign  of  tangential 
components  may  be  used  for  any  second-order  symmetric  tensors.  It  is  a complete  mathematk 
tool  that  is  compatible  with  the  analytical  method,  as  shown  in  the  following  applications. 
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Figure  2.  Stress  characteristic  directions. 


A convention  is  made  in  (2)  that  the  upper  and  lower  of  the  double  sign  T or  ± are  related  to  the  +m 
and  -m  characteristic  directions,  respectively.  This  convention  is  observed  throughout  the  paper. 

Let  the  a<oordinate  of  the  center  C and  the  radius  of  the  Mohr  circle  be  o and  (a),  respectively. 
Then,  the  stress  components  are  given  by 

= o + @ cos  20 

Oy  = o - © cos  26  (3) 

T„y  = © sin  26. 

The  yield  criterion  is  given  by 

© = asinp  + ccosp  (4) 

where  c is  the  cohesion.  ® is  always  positive,  but  a can  be  negative,  zero,  or  positive. 

Equation  (1)  may  be  analytically  derived.  This  can  be  achieved  by  substituting  (3)  and  (4)  into 
the  stress  balance  equations  and  applying  the  theory  of  characteristic  lines  (Courant  and  Hilbert,^ 
Abbott' ).  This  method  is  illustrated  later  by  the  derivation  of  (18).  Sokolovski'*  shows  a slightly 
different  but  useful  approach  in  deriving  (1 ).  (Although  Sokolovski’s'*  central  angle,  which  he 
denotes  by  20,  is  < BCV  in  terms  of  Figure  2,  his  formulas  [ref.  18,  (1.15)  and  (1.16)]  are  the  same 
as  (1)  in  this  paper.) 


ANALYSIS  OF  STRAIN  RATE 
Constitutive  equations 

We  shall  discover  the  constitutive  equations  of  the  Coulomb  material  by  combined  use  of  the 
theory  of  characteristic  lines  and  the  geometry  of  the  Mohr  circle. 
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We  consider  a strain-rate  Mohr  circle  determined  in  response  to  the  stress  Mohr  circle.  Let  us 
superimpose,  as  shown  in  Figure  3,  the  strain-rate  Mohr  circle  on  the  stress  Mohr  circle,  by,  if 
necessary,  changing  the  scale,  moving  the  origin,  and  rotating  the  axis  of  the  strain-rate  Mohr  circle, 
so  that  the  pole  of  the  latter  falls  on  the  pole  of  the  former.  In  Figure  3,  the  axis  A^CB^  of  the 
strain-rate  Mohr  circle  is  rotated  by  angle  2x  from  the  axis  A„  CB,,  of  the  stress  Mohr  circle,  where 
X denotes  the  non-coaxial  angle.  The  origin  of  the  e,  y coordinates  is  at  0^.  The  strain-rate  Mohr 
circle  of  zero-length  radius  is  excluded  from  this  operation. 

First,  let  us  formulate  the  strain  rate_components.  Let  the  6<oordinate  of  the  wnter  C and  the 
radius  of  the  strain  rate  Mohr  circle  be  e and  respectively;  (D  is  positive,  but  e can  be  positive, 
zero,  or  negative.  The  values  of  and  are  given  by  the  intersection  of  axis  A^  CB^  with  the 
normals  drawn  from  V and  H to  A^  CB^,  respectively.  Then,  by  using  the  relation  < VCA^  = 26, 
we  find 

®xx  = f © coi(26  + 2x) 

^yy  = f - © cos(2e  + 2x)  (5) 

«xy  = © + 2x). 

Let  and  be  the  x-  andy-  components  of  the  velocity  v.  We  define  the  compression  to  be 
positive;  'Hen, 

• = ^ 

~ dx 
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Second,  we  introduce  a lemma:  The  necessary  and  sufficient  condition  for  the  x-direction  to  be 
an  4<haracteristic  direction  is 

c,  = 0.  (7) 

To  prove  this,  we  assume  the  constitutive  equations  in  the  most  general  linetr  form 

^ii  = + (8) 

by  applying  the  summation  rule  for  the  repeated  indexes,  h and  k,  where  e’jj  represents  one  of  the 
components  e^y,  or  e^y,-  Ohk  represents  one  of  the  components  o^,  Oy  or  r,,y;/4**  and  fijj  are 
constants;  and  X is  a nonzero  constant. 

Let  us  replace  the  left-hand  side  of  (8)  with  the  right-hand  side  of  (6),  and  consider  the  equations 
thus  found  to  be  the  differential  equations  for  determination  of  and  Vy,  assuming  that  all  the 
quantities,  except  X,  on  the  right-hand  sides  of  (8)  are  given.  Let  us  suppose  that  the  x-direction  is 
the  characteristic  direction  of  the  system  of  differential  equations  in  (8);  i.e.,  we  suppose  that  the 
x-direction  is  such  that,  given  and  Vy  along  the  x-direction,  we  cannot  determine  dvjby  or 
dv^jdy  by  substituting  the  known  functions  into  the  differential  equations  in  (8)  (Courant  and 
Hilbert,’  Abbott'). 

If  is  not  zero,  the  equation  of  e'^x  in  (8)  can  be  used  to  determine  X,  because  X is  the  only 
one  unknown  in  this  equation.  If  (7)  is  correct,  the  right-hand  side  of  the  equation  of  in  (8) 
must  be  equal  to  zero,  because  this  equation  is  an  identity;  therefore,  X is  indeterminate.  Then, 
bvjby  and  bv^ldy  are  also  indeterminate.  Equation  (7)  is  therefore  the  necessary  and  sufficient 
condition.  The  lemma  is  thus  proved.  We  shall  later  show  that  we  can  actually  arrive  at  the  constitu- 
tive equations  in  the  form  of  (8). 

Third,  we  intioduce  the  principie  of  partiai  coincidence  as  a requirement  for  a plastic  deforma- 
tion to  occur.  In  the  conventional  theory  of  plasticity,  it  is  postulated  that  both  of  the  stress- 
characteristic  directions  are  coincident  with  both  of  the  strain-rate  characteristic  directions.  How- 
ever, the  coincidence  of  both  of  the  characteristic  directions  is  not  necessary.  We  postulate  that, 
for  a piastic  deformation  to  occur,  one  member  of  a set  of  directions  must  be  coincident  with  the 
member  of  the  same  sign  m in  the  other  set.  If  one  member  is  coincident,  the  boundary  flow  condi- 
tions and  the  boundary  stress  conditions  can  be  met  along  the  coincident  direction.  If  there  is  no 
coincident  direction,  it  is  impossible  to  connect  the  two  sets  of  discontinuous  solutions  across  a line. 

The  coincident  and  noncoincident  directions  are  called  doubiet  and  singlet,  respectively.  There 
are  two  singlets;  a o-singlet  and  an  e-singlet.  When  the  doublet  has  sign  +m  or  ~m,  we  say  that  the 
+m  or  -m  coincidence  has  occurred.  Finally,  if  velocity  is  in  the  doublet  direction,  we  call  the  flow 
line  a slip  line. 

Fourth,  let  us  resume  the  interpretation  of  (7).  We  assume  that  the  x-direction  is  the  doublet; 
then,  in  accordance  with  the  +m  or  -m  coincidence,  the  x-direction  must  be  either  PD  or  PE  in 
Figure  2.  Then,  from  the  geometry,  6 must  be  equal  to  ±6.  Substituting  (7)  and  the  above-found 
value  of  0 into  (5) , , we  find 


e = - ©cos(2x  ± 26). 


(9) 


Therefore,  (9)  is  one  of  the  necessary  conditions  for  a plastic  deformation  to  occur. 

We  shall  show  that  we  have  now  a set  of  sufficient  conditions  for  a plastic  deformation  to  occur. 
Equations  (3),  (4),  (5)  and  (9)  may  be  regarded  as  a parametric  expression  of  the  constitutive  equa- 
tions. Eliminating  parameters  6,  o,  and  e,  we  can  find  explicit  forms  of  the  constitutive  equations. 
The  simplest  is  the  following  matrix  expression; 


(10) 
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where  a prime  si^ifies  that  the  respective  quantity  is  deviatoric.  This  equation  shows  that  X in 
(8)  is  given  by  ©/©.  If  the  deviators  C”  = Cjj-e  and  o|j  = Ojj  - o are  relegated  to  restore  the 
regular  components  e'lj  and  a|j,  we  find,  as  shown  in  ref.  20,  that  (10)  yields  the  most  general  form 
(8)  Note  that  the  constitutive  equations  thus  derived  contain  x- 

Equation  (9)  is  as  significant  to  the  strain-rate  tensor  as  the  Coulomb  criterion  (4)  is  to  the  stress 
tensor,  and  is  called  the  volume  characteristic.  The  rate  of  volurrtc  change 


e 


2 \ 9x 


(11) 


can  be  positive,  zero,  or  negative  depending  on  the  value  of  non-coaxial  angle  x.  The  volume 
characteristic  was  derived  in  ref.  20  without  assuming  (8),  but  the  use  of  the  theory  of  characteristic 
lines  in  this  paper  is  much  more  elemental  than  in  ref.  20.  It  was  shown  in  ref.  20  that  the  theory 
of  plastic  potential  does  or  does  not  apply  when  either  x = 0 or  # 0,  respectively. 

Density  y can  be  determined  by  use  of  i : the  equation  of  the  conservation  of  mass 


+ div(7v)  = 0 
or 


yields 

la  + V.  =2ly 
dt  ds 


(12) 


where  v is  the  velocity  vector,  v the  velocity  magnitude,  and  s the  length  of  the  flow  line. 
Strain-rate  characteristic  directions 

We  shall  write  out  the  e-singlet  direction  which  is  still  not  formulated.  Eliminating  the  two 
parameters  e and  © from  the  three  equations  in  (5),  we  find 


A 


dv^  dv^  dv^  di/„ 


= 0 


(13) 


where 


A = sin  (20  + 2x) 
B = cos(20  + 2x). 
We  rewrite  (9)  to 


bv 


-J1  + r iik 
bx  ^ by 


+ C 


bVy 

bx 


+ A 


bVy 

by 


= 0 


(14) 


where 


C = cos(2x  i 26) 
by  replacing  e with  (14)  and 


. av  , ^v 

dv  = dx  + 

* Ox  by 


eliminating  © by  use  of  (5)3.  Two  uniqueness  conditions 
dy 


(15) 
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3y 


(16) 


dx  + 


must  also  be  satisfied. 

Regard  equations  (13),  (14),  (1 5),  and  (16)  as  the  simultaneous  linear  equations  with  four  un- 
knowns; bvjbx,  bv^jbx,  dvjby,  and  bv^jdy.  We  can  determine  the  characteristic  directions  by 
letting  the  denominator  determinant  equal  zero  (Courant  and  Hilbert,^  Abbott' ) 


A _B  -B  -A 
A C C A 
dx  dy  0 0 

0 Q dx  dy 


= 0. 


(17) 


The  quadratic  equation  (17)  has  two  roots:  one  is  (2),  the  doublet;  the  other 's  the  e’-singlet 

^ = tan(0  + 2x±6).  (18) 

Equations  (17)  is  the  denominator  of  the  solutions  of  the  simultaneous  linear  equations  (13), 

(14),  (15),  and  (16)  for  the  unknowns  bvjbx,  dv^ldy,  dv^jbx,  and  dVy/dy.  Then,  the  numerators 
of  the  rational  expressions  of  these  unknowns  must  also  be  equal  to  zero.  We  can  show  that  all  four 
numerators  reduce  to  a single  equation 

dv^dx  + dVy  dy  = 0 (19) 

where  dyjdx  is  either  (2)  or  (18). 

If  we  simply  want  to  verify  (19),  it  is  much  simpler  to  form  a linear  combination  of  (13)  and  (14). 
We  can  find  that  both  (1 3)  x cos(20  + 25)  + (14)  and  (1 3)  x cos(20  + 4x  ± 26)  + (14)  reduce  to  the 
same  formula 

dVv  , \ 3i/„  , 

ST  ^ 1^'j  ^ ° 

where  dyidx  is  either  (2)  or  (18).  Equation  (20)  reduces  to  (19)  by  use  of  (15)  and  (16).  Equation 
(20)  is  equivalent  to 

«’xx  2exy  dx  dy  + (dy)'^  = 0.  (21 ) 

Therefore  (19)  means  that  the  strain-rate  characteristic  lines  neither  elongate  nor  shrink. 

Constitutive  geometry 

It  is  possible  to  express  the  constitutive  relationships  with  the  geometry  of  the  +m  and  -m 
coincidences,  as  shown  in  Figures  4a  and  4b,  respectively. 

In  Figures  4a  and  4b,  and  -m  doublets  are  at  D and  E,  +m  and  -m  a-singlets  are  at  E^,  and  D^, 
and  +m  and  -m  e-singlets  are  at  E^  and  D,,  respectively.  Equations  (2)  and  (1 8)  give,  as  shown  at 
P,  the  angles  that  a- and  e'-characteristic  directions  make  with  the  x-axis.  Point  0^  is  the  intersection 
of  the  two  tangents  drawn  at  the  strain-rate  characteristic  points.  Details  of  the  geometric  relation- 
ship are  described  in  the  following. 

Proposition  1.  The  e-singlet  direction  makes  angle  2x  ± 28  with  the  doublet  direction. 

Proof.  Compare  the  arguments  of  the  tangent  functions  in  (18)  and  (2). 

Proposition  2.  The  e'-singlet  is  by  the  central  angle  4x  apart  from  the  a-singlet. 

Proof  In  the  case  of  Figure  4a,  note  that 
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Proposition  3.  The  e-singlet  is  symmetric  with  the  doublet  with  regard  to  the  axis  A,  B,. 

Proof.  !n  the  case  of  Figure  4a,  note  that 

< Dco;  = < Dco;  + < 0;c0;  = 26  + 2x 

< o;CE,  = -<  0;C0;  + < 0;CE,  + < E„CE,  = 26  + 2X. 

Therefore  D and  E,  are  symmetric  with  regard  to  B^CA, . In  the  case  of  Figure  4b,  note  that 

< O^CE  = < 0;CE-<  0;c0;  = 26 -2X 

< D,CO;  = < D„CO;-<  D„CD,  + < 0;c0;  = 26 -2X. 

Therefore  E and  are  symmetric  with  regard  to  B^CA,. 

Proposition  4.  Point  is  on  the  line  A^CB^. 

Proof.  This  is  a corollary  of  Proposition  3, 

Proposition  5.  The  strain  rate  counterpart  of  the  frictional  angle  is  p - 2x  in  the  case  of  Figure 
4a  and  p + 2x  in  the  case  of  Figure  4b. 

Proof.  The  strain  rate  counterpart  of  the  frictional  angle  is  < DO^C  in  the  case  of  Figure  4a  and 
< DjOjC  in  the  case  of  Figure  4b.  They  are  given  by  p - 2x  and  p + 2x,  respectively. 

Proposition  6.  The  origin  of  the  e,  7 coordinates  is  at  the  point  at  which  A^  B,  intersects  with 
DEj  in  the  case  of  Figure  4a  and  with  E in  the  case  of  Figure  4b.  The  direction  of  the  e-axis  must 
be  determined  to  conform  with  (9). 

Proof.  We  can  prove  in  the  case  of  Figure  4a 

< Dco;  = < Dco;  + < o;co;  = 2x  + 26 

and  in  the  case  of  Figure  4b 

< Eco;  = < ECO;-<  o;CO;  = 26-2x. 

Therefore,  use  of  (9)  shows  that  must  be  located  as  stated  in  the  proposition.  The  directions  of 
the  e‘-axes  shown  in  Figures  4a  and  4b  are  chosen  to  conform  with  (9)  for  the  values  of  2x  as  shown 
in  these  figures.  Note  that  the  negative  and  positive  values  of  e are  defined  in  this  paper  as  dilation 
and  contraction,  respectively. 

Proposition  7.  In  the  case  of  +m  coincidence,  e is  positive  or  negative  when  A^  in  Figure  5a  is  on 
the  arcs  MDN  or  ME„  N,  respectively.  In  the  case  of  -m  coincidence,  e is  positive  or  negative  when 
A^  in  Figure  5b  is  on  the  arcs  MEN  or  MD„N,  respectively.  Line  MCN  in  Figure  5a  or  5b  is  drawn 
parallel  to  Oj,  D or  E,  respectively. 

Proof.  Use  of  propositions  6 and  3 shows  that,  when  A is  at  M,  C is  the  origin  of  the  i,  y co- 
ordinates. On  which  side  of  MN  point  A^  must  be  located  to  make  e'  > 0 or  < 0 can  be  found  by 
use  of  proposition  6 or  (9). 


STRAIN-RATE  TENSOR 
The  dyadic  expression 

The  relationships  of  the  strain  rate  components  in  (5)  and  (9)  give,  as  shown  below,  a simplified 
expression  of  strain  rate  tensor,  which  facilitates  the  transformation  of  curvilinear  coordinates.  To 
show  this,  second-order  tensors  must  be  expressed  as  dyadics  (Wilson,”  Brand,’  Sedov,“  Yoshimura,** 
Takagi”). 

Letting  the  x-  and  y-directions  be  c^and  Cy , respectively,  the  strain  rate  tensor^  may  be  expressed 
as 
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Figure  5a.  Domains  e % 0 on  the  Mohr 
circle.  The  case  of  Hn  coincidence  Is 
shown. 


^ ^XX  ®X  ^X  ®XV  ^^X  ^'V  ^'V  ®x)  ^'v 


Figure  5b.  Domains  e ^ Oon  the  Mohr 
circle.  The  case  of  -m  coincidence  is 
shown. 


C„ 


(22) 


The  bases  c^,  Cy,  Cy  Cy,  and  Cy  Cy  of  this  tensor  are  dyads  (Wilson,*’  Brand*).  (Note  that 
Cx  Cy  * Cy  i.e.,  vectors  are  not  commutative  in  a dyadic  product.)  In  the  conventional  tensor 
analysis,  tensor  bases  being  omitted,  the  components  are  transformed;  however,  use  of  invariants  e 
and  © forces  us  to  transform  tensor  bases.  Substituting  (5),  we  transform  (22)  to 

^(Cx'x  + S Cy)  + 

© j [Cx  cos(0  T 6)  + Cy  sin(0  T 5)]  [c,<  cos(0  + 2x  ± 5)  + Cy  sin(e  + 2x  ± 5)1  - 

- (c„  sin(0-T  6)  - Cy  cos(0  ? 6)]  [c^  sin(0  + 2x  ± 6)  - Cy  cos(0  + 2x  ± 6)i  j . (23) 

This  will  further  be  transformed  in  the  following  to  a compact  dyadic  expression. 

Case  1.  When  2x^28  * nn  (n  is  an  integer),  the  doublet  (2)  and  the  e-singlet  (18)  are  distinct. 
We  denote  the  directions  of  the  doublet  and  the  e-singlet  by  the  subscripted  unit  base  vectors  u, 
and  U2,  respectively. 


u,  = c^  cos(6  T 5)  + Cy  sin(0  + 5) 

“2  ~ ‘^x  ^o*(®  + 2x  ± 5)  + Cy  sin(0  + 2x  ± 5). 


(24) 


When  2x  ± 26  # 0,  vectors  u^  and  U2  are  not  orthogonal;  then  it  is  convenient  to  introduce 
superscripted  vectors  e'  and  e^  defined  by  ti-  Uj  = 6j,  where  6|  is  a Kronecker  delu,  and  a-  b 
means  the  dot  (scalar)  product  of  vectors  a and  b.  It  is  more  convenient  in  this  case  to  use  the 
directions  u'  and  u^  of  e'  and  e^,  respectively,  instead  of  e'  and  e^  themselves;  in  more  detail, 
components,  if  appropriately  chosen,  of  a tensor  expressed  with  dual  unit  vectors  U|,  u'  (/  = 1 , 2) 
can  mean  physical  components,  i.e.,  components  with  clear  physical  meaning.  We  determine  the 
senses,  which  are  still  left  undecided,  of  the  superscripted  unit  base  vectors  u' , u^  by  making  them 
satisfy 


u,u’  = U2u2 


Thus,  the  two  vectors  arc  determined: 

u’  = c,<  sin(0  + 2x  ± 6)  - Cy  cos(0  + 2x  ± 6) 
u^  = - sin(0  + 6)  + Cy  cos(0  T 6). 


(25) 


(26) 
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Us«  of  (24)  and  (26)  transforn's  (23)  to 


i 


e = e(CxCx+c^Cy)+  0(u,  U2  + u2u').  (27) 

Solving  (26)  for  and  c^,  we  get 

sin(2x  i 26)  = u’  co$(0  ? 6)  ♦ cos(0  ♦ 2x  ± 6) 

(28) 

$in(2x  i 25)  = u'  sin(0  + 6)  + sin(0  + 2x  ± 25). 

Substituting  (28)  into  (24),  we  get 

u,  sin(2x  ± 26)  = u'  + cos(2x  ± 26) 

Uj  sin(2x  ± 26)  = u'  cos(2x  ± 26)  + u^. 

By  use  of  the  formulas  derived  above,  we  can  transform  all  the  base  vectors  in  6 in  (27)  to  u'  and  u^; 
then,  £ should  take  the  form 

£ = e,,  u’  u’  + e,2  (u’  u' ) + £22  • (29) 

Actual  calculation  yields  that 


e,,  = 0 

€22  = 0 (30) 

ei2  = ©■ 

The  actual  calculation  shows  that  (30) ^ and  (30)2  are  the  two  equations  in  (19),  as  shown  by  (21). 

Vector  u'  defined  by  the  convention  (25)  does  not  necessarily  make  an  acute  angle  with  Uj,  as  in 
the  case  of  the  orthodox  definition  e'.ej  = 6j.  However,  the  angles  between  Uj  and  o'  are  simul- 
taneously either  acute  or  obtuse;  then,  dyads  in  (29)  do  not  change  sign;  therefore,  the  components 
do  not  change  sign.  Thus,  the  convention  (25)  completely  serves  our  purpose.  Finally,  let  us  note 
that  we  prefer  the  new  terminology  of  "subscripted"  and  "superscripted”  to  the  classical  terminology 
of  "covariant”  and  "contravariant.” 

Case  2.  When  2x  ± 26  = rttr,  there  is  a single  coincident  6-characteristic  direction;  and  the  de- 
composition (23)  does  not  work.  In  this  case,  (9)  and  (5)  reduce  to 

6 = (-1)"-i  © 

^xx  " (-I)""’ 2©sin7  (0+6) 
eyy  = (-!)"■'  2 ©cos7  (6  + 5) 

^xy  = (-1)”  ©sin(26  + 26). 

Then,  (22)  transforms  to 

£ = (-I)"’’  2 © [c„  sin(6  + 6)  -Cy  cos(6  ? 6)1  |c„  sin(6  ? 6)  - Cy  cos(0  T 6)). 

We  will  transform  this  to  a compact  dyadic  expression. 
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Let  Uf  and  Uj  be  the  directions  of  the  doublet  and  the  o-singlet,  respectively, 


r 


! 


'i 


Uf  = c„  co$(fl  T 6)  + Cy  sin(0  T 6) 

Uj  = Cx  co$(0  ±6)  + Cy  sinjfl  ±6). 

Let  and  be  the  reciprocal  directions  | 

u’  = Cj,  sin((i  T6)  - Cy  cos(0  ±6)  | 

= -c^  sin(6  T 6)  + Cy  co5(fl  T 6)  ] 

satisfying  the  convention  in  (25).  Then  (31)  becomes  | 

e=  (-1)"-’ 2 0u2u2.  (32) 

Therefore,  components  of  ^ are  given  by 
*11  = 0 

622  = (-1)"-’2  © (33) 

^12  = 0- 

When2x±26  = nw,  (1 9)  yields  only  one  equation,  i = 0.  The  missing  equation  is  supplied  here; 

itis€j2  = 0. 

Plastic  work  rate 

The  radius  © is  a factor  of  the  plastic  work  rate.  The  plastic  work  rate  W is  defined  by 

\V  - *xx  ^^xy  *xy  *^y  ^yy  • 

Substituting  (3)  and  (S)  and  using  (4)  and  (9),  it  becomes 

iV  * 2 © cot  p [c  cos(2x  ± 26)  ± © sin  2x1 . (34) 

This  quantity  must  be  positive  or  zero,  if  no  other  irreversible  process  is  concurrent.  The  condition 
© i 0 (35) 

must  always  be  satisfied. 

Coordinate  transformation 

Equations  of  motion  are  found  by  transforming  either  (19)  or  (33)3.  In  addition  we  must 
formulate  © to  evaluate  density  7 by  use  of  (1 2).  In  most  cases,  the  coordinates  of  (19)  can  easily 
be  transformed.  However,  formulation  of  €^2  (33)3,  and  ^22  (33)2. 

velocity  components  needs  elaborate  transformation  of  the  strain-rate  tensor,  as  described  below. 

Because,  in  the  following  we  must  deal  with  the  skew  curvilinear  coordinate  system,  we  need  a 
more  expressive  tensor  notation  than  the  conventional.  We  express  the  strain  rate  tensor^  in  the  fol- 
lowing invariant  differential  form.  Letting x,y,  c,,,  and  Cy  be  renamed  f’,  e^,  and  e^,  we  rewrite 

e,  defined  by  (22),  by  use  of  (6),  as 
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e ^ .1 L ^ ^ ^ t' 
2 [ af'  3^’ 


where  / is  the  summation  index  and 
‘^x  ■*■  ‘'v  S • 

The  operator 

ei-^  = c,  i-+c,  f (38) 

g|i  * dx  r dy 

is  an  invariant,,  denoted  by  V,  or,  as  often  called  grad.  Equation  (36)  is  convenient  for  complicated 
coordinate  transformation. 

We  define  curvilinear  coordinates  and  by 

X = x(?',^2) 

(39) 

y y(l',?2). 

Curves  I ' and  are  such  curves  on  which  only  f ' and  respectively,  vary  and  and  f’  are, 
respectively,  kept  constant.  We  define  e,  and  e2  by  vectors  tangent  to  curves  and  respectively. 
Lengths  of  ej  and  62  must  satisfy  the  relation 

dr  = c^dx  + Cy  dy 

(40) 

= e,  dV  +e2d|2. 

This  equation  gives  e,  and  62  as  follows; 

e,  = JL  = c,  iif.  (41) 

ai'  3^'  9?‘ 

We  then  define  reciprocal  vectors  e'  by 

e'.e|  = 6j.  (42) 

We  may  reinterpret  the  operator  e'  3/31'  on  the  left-hand  side  of  (38)  in  terms  of  the  new  defini- 
tion. We  can  prove  that  (38)  is  still  true  even  in  the  new  definition.  To  prove  this,  note  that  e^  may 
be  written  as 

because  this  expression  satisfies  (42).  Substitution  of  (43)  reduces  the  left-hand  side  of  (38)  to  the 
right-hand  side. 

Note  that  the  derivaties  3e|/3|*  and  3e*/3f'  can  be  calculated  by  use  of  (41 ) and  (43),  respectively, 
because  c„  and  Cy  are  constant.  Therefore,  even  if  v may  be  expressed  with  bases  e,-  or  e*,  the  right- 
hand  side  of  (36)  can  be  reduced  to  the  right-hand  side  of  (29),  where  Uj  and  u*  are  directions  of 
Cj  and  e',  respectively. 


Figure  56,  Geometry  of  the  chavcteristic  directions  of  Figure  6a. 
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EXAMPLE  ' 

J 

Figure  6a  shows  a net  of  stress  characteristic  lines  caused  in  the  ground  sustaining  a vertical 
rectangular  load  applied  along  the  segment  da  on  the  horizontal  ground  surface.  The  line  abpcd  is  ! 

the  slip  line.  Hill"  gave  the  original  solution,  as  explained  by  Prager  and  Hodge,"  for  the  case  of 
frictional  angle  p = 0.  Figure  6a  is  an  extension  of  Hill's  solution  to  the  case  of  p # 0,  as  explained  i 

in  the  following.  We  shall  fit  the  simplest  velocity  solution  into  the  extended  Hill's  solution. 

The  stress  solution 

In  terms  of  the  geometry  of  the  Mohr  circle,  the  solution  may  be  explained  as  follows:  ' 

In  Figure  6a,  the  soil  below  the  line  abpcd  does  not  move  and  only  the  soil  above  abpcd  moves.  i 

We  may  assume  that  the  tangential  stress  on  the  upper  side  of  the  line  abpcd  is  in  the  direction  of  the  j 

soil  flow;  then,  the  sign  m on  the  line  abpcd  is  + . Therefore,  the  doublet  is  at  point  D in  Figure  6b,  j 

and  the  O'Singlet  is  at  . 
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In  the  triangular  regions  oab  and  ocd  (Fig.  6a),  the  stress-characteristic  lines  are  rectilinear.  In 
the  former,  the  pole  is  at  point  B (Fig.  6b);  the  doublet  direction  BD  and  the  o-singlet  direction 
BE„  are  parallel  to  ba  and  bo,  respectively.  In  the  latter,  the  pole  is  at  point  A;  the  doublet  direction 
AD  and  the  o-singlet  direction  AE„  are  parallel  to  cd  and  co,  respectively.  The  former  is  called  the 
active  Rankine  sute,  the  latter  the  passive  Rankine  state  (Terzaghi**). 

In  the  region  obc,  the  pole  P moves  from  B to  A through  the  upper  semicircle,  as  point  p moves 
from  b to  c.  The  tangent  at  point  p is  parallel  to  PD,  and  the  radial  direction  op  is  parallel  to  PE„. 

The  o-singlets  are  a pencil  of  straight  lines  passing  through  point  o.  

Let  ^ be  the  angle  designating  the  position  of  a o-singlet  op  measured  from  ob  in  Figure  6a 

< bop  = <1/ . 

Then,  in  Figure  6b,  we  find 

< BE„P  = 

Thus,  we  can  prove  that 


0 + i)/  = 7r/2 


(E-1) 


because,  as  shown  in  Figure  6b,  < PE„  A = 0.  Eliminating  0 in  (2)  by  use  of  (E-1 ),  we  get  the 
equation  of  the  o-singlet  lines  in  Figure  6a, 

y = X cot(i)/  - 6).  (E-2) 

We  can  otherwise  derive  (E-2)  directly  from  the  geometry  of  Figure  6a,  because  < boy'  = 6, 
where  ^ is  the  downward  extension  of  the  y-axis. 

The  doublet  is  given  by 

dyjdx  = cot(\I/  + 6}  (E-3) 

We  can  transform  this,  by  eliminating  \j/  by  use  of  (E-2),  to  a total  differential  equation 
{x  dy  - y dx)  + (x  dx  +>'  dy)  tan  26  =0. 


To  integrate  this,  let 

X = - r sin(\|/  - 6) 
y = -rcos(^  -6) 


(E-4) 


where  r is  the  radial  coordinate  of  a point  under  consideration.  Then  we  get 

(E.5) 

where  f is  the  value  of  r at  ^ = 0,  i.e.,  on  the  initial  line  ob. 

Velocity  equations  in  the  o-characteristic  curvilinear  coordinates 

We  shall  fit  a velocity  field  into  the  stress  field  obtained  above.  Because  x unknown,  we  shall 
use  the  o-characteristic  curvilinear  coordinates  Because  this  curvilinear  coordinate  system  is 
skew,  we  need  a more  general  and  simpler  approach,  i.e.,  use  of  base  vectors,  to  Introduce  skew 
curvilinear  analysis  Let  the  directions  of  curves  and  be  u^  and  u^. 
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Uf  = -c^  sin(^  - 6)  - Cy  co$(^ -6)  (E-6) 

= - Cx  sin(^  + 8) -Cy  cos(^  + 6).  (E-7) 

These  equations  are  found  by  use  of  (2)  with  the  substitution  of  0 from  (E-1). 

In  addition  to  the  subscripted  base  vectors  Uj  and  u^,  we  introduce  the  superscripted  base  vectors 
uf  and  u* . They  are  unit  vectors  defined  by 

= 0 

uf.uj.  = 0 

Uj. u^  = u^.u*  = sin  26. 

They  are  given  by 

uf  = c^  cos(^  + 6)  - Cy  sin(^  + 6)  (E-8) 

u'l'  = -c^  cos(^ -6)  + Cy  sin(^ -6).  (E-9) 

These  four  vectors  are  directed  as  shown  in  Figure  6a.  The  radial  op  makes  angle  26  with  the  doublet 
at  p,  because  Uj..u^  = cos  26. 

We  decompose  velocity  v into  the  directions  u^  and  u^, 

V = V,  u^  + V2  (E-10) 

Note  that  u?  is  normal  to  u^ . Substituting  from  (£-7)  and  u^  from  (E-8)  into  (E-10),  and  com- 
paring the  result  with 

= '^x  Cx  ‘'y  'v 


‘'x  ~ “ ‘'l  sin(V'  + 6)  + ^2  cos(^  +6) 

Vy  = - V,  cos(4i  + 6)  - ^2  sin(^  + 6). 

Then,  substituting  and  Vy  from  (E-1 1 ) and  dyidx  from  (E-3),  (1 9)  in  the  case  of  the  doublet 
direction  transforms  to 


(E-IJ 

where  s is  the  length  along  a doublet  curve.  The  e-singlet  is  +m  coincident  and  is  given  from  (1 8)  by 

dyjdx  = cot(^-2x-6).  (E-1 3 

Substituting  and  Vy  from  (E-1 1 ) and  dyjdx  from  (E-1 3),  (1 9)  in  the  case  of  the  e’-singlet  direction 
becomes 


cos(2x  + 26)^^  ^ ‘'2  |j^)  - «n(2x  + 26)^^  - v,  = 0 


where  t is  the  length  along  the  e'-singlet  curve. 
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We  shall  express  (E-12)  and  (E-14)  with  curvilinear  coordinates  \1j,  f.  Coordinates  \l/  and  f are 


related  to  coordinates  x,  y and  s,  t by 

dr  = c„  cfx  + Cy  dy 

(E-1 5) 

+ 

11 

(E-16) 

= u^ds  + u,  dt 

(E-1 7) 

where  and  u,  are  directions  of  the  doublet  and  the  e-singlet,  respectively.  We  cannot  use  the 
e-characteristic  curvilinear  coordinates  in  s-  and  f-dircctions,  because  x is  unknown. 

We  express  x,  y as  functions  of  f and  \l/  by  substituting  (E-5)  into  (E-4),  substitute  these  into 
(E-1 5),  and  compare  the  coefficients  of  df  and  d>j/  in  the  transformed  (E-1 5)  and  the  untransformed 
{E-16);  thus  we  find 


ej  = e'^' 


p = — i lanp  y 

sin  26 


(E-1 8) 


It  is  obvious  that 

Uj  = Uf.  (E-1 9) 

In  the  Cartesian  coordinates,  u,  is  given  by 

u,  = - c,,  sin{t;/  - 2x  - 6)  - Cy  cos(i//  - 2x  - 5)  (E-20) 

where  (E-1 3)  is  used.  We  solve  (E-6)  and  {E-7)  for  c^  and  c^,  substitute  these  results  into  (E-20),  and 
find 


Uj  sin  26  = Uj.  sin(2x  + 26)  - sin  2x.  (E-21) 

We  substitute  (E-1 8)  into  (E-1 6),  and  (E-1 9)  and  (E-21 ) into  (E-1 7),  compare  the  coefficients  of 
Uj.  and  u,^  in  the  two  transformed  equations,  and  find 

gi^  tan p = sin(2x  26) 
sin  26 


tin  p dt. 

sin  26  sin  26 


Use  of  the  above  equations  yields 


pii  tanp  ^ = sin  2o  3 
3s  5-  3i)/ 


pi;/ tanp  A - sin(2x  -I-  26)  3 sin  2x  3 
3f  sin  26  3r  " f 3i/t  ‘ 


Thus  (E-12)  becomes 
3vi 


'1 


dtp 


+ V/-  = 0 


(E-22) 


(E-23) 
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r 


and  (E-14)  becomes 

- 5755  ^ ■ “■ 


(E-24) 


The  constant  speed  solution 

We  shall  find  the  simplest  velocity  field  that  can  ^fitted  in  Figure  6a.  If  we  denote  by  V the 
velocity  on  the  segment  Sa  in  the  doublet  direction  ab,  the  simplest  solution  in  the  rectilinear  region 
oab  under  the  condition  that  oab  is  the  slip  line  is 


*'1 

‘'2 


0 


(E-25) 


where  the  decomposition  of  v into  the  doublet  direction  and  the  normal  to  the  doublet  direction  as 
shown  in  (E-10)  is  still  used.  Equations  (E-25)  are  the  solution  even  in  the  curvilinear  region  obc,  if 


X = 0 


(E-26) 


because  (E-2S)  satisfies  (E-23)  and  (E-24)  on  this  assumption.  We  can  find,  if  we  will,  more  com- 
plicated solutions;  however,  in  this  paper,  we  are  not  interested  in  deriving  them. 

In  the  rectilinear  regions  oab,  the  above  solution  gives  (j)  = 0,  and  therefore  W = 0.  However, 
this  conclusion  does  not  hold  true  in  the  curvilinear  region  obc.  To  calculate  0 in  this  case,  we 
shall  compute  e'  (3v/a|')  contained  in  (36),  which,  in  this  case,  becomes 


eiiv-  = ef|f + e*  IV- 
d|i  3'/' 


(E-27) 


The  superscripted  base  vectors  are  determined  by  applying  (E-1 8)  to  (42) 


ef  = _L_  e-*  uf 

cosp 


= L e"^  u'^  . 


(E-28) 


To  differentiate  v in  (E-10),  note  that  (E-7)  and  (E-8)  yield  the  following  relations: 


^ = 


-uf 


di// 

lf-» 


3uf  _ u 


Thus  (E-27)  transforms  to 
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g\j/  tanp 


cosp  e 


i Jv 

a^i 


(E-29) 


where  we  have  used  (E-28)  to  change  e^  and  e'*'  to  and  u'^,  respectively. 

We  shall  transform  (E-29)  to  the  form  of  (29);  i.e.,  we  shall  change  the  base  vectors  in  (E'29)  to 
u*  and  u',  where  Uj  and  u,  are  given  by  (E-19)  and  (E-20),  respectively.  We  find 


= - c^  cos(i;'  - 2x  - 5)  + Cy  sin(V/  - 2x  - S) 
u'  = c^  cos(\p  + 6)  - Cy  sin(^  + 5) 


where  the  condition  (25)  is  observed.  Solving  the  above  two  equations  for  and  Cy,  and  substi- 
tuting this  solution  into  (E-7),  (E-8),  and  (E-9),  we  can  express  u^,  uf,  and  u'*'  in  terms  of  u*  and 
u'.  Substituting  these  results  into  (E-29),  we  find 


e^  tan  p cosp 


sin^(2x  + 26) 


= In 


u'  sin(2x  + 26)  [u*  + u'  cos(2x  + 26)] 


+ 


+ u'  u^  sin^(2x  + 26)  + 

+ A!E2S  ^ s'f  26  - u^  sin  2xj  fu^  + u^  cos(2x  + 26)J  + 

+ - V,  j [u^  sin  26  - u*  sin  2x1  u'  sin(2x  + 26) . 


(E-30) 


Changing  the  order  of  base  vectors  in  the  dyadic  products  in  (E-30),  we  construct  the  right-hand 
side  of  (36);  then,  we  can  determine  coefficients  €,,,622,  and  6,2  in  (29).  Thus,  from  the  expression 
of  e,2,  we  find 


(d)  = i e't^  “np  cosec(2x+  26) 


sin^  26 /.,  ^‘'2\ 


The  expressions  of  the  components  e,,  and  €'22  give  equations  (E-23)  and  (E-24),  respectively. 
Therefore,  in  the  simplest  solution  (E-25),  © is  given  by 


© = i 1/  cos  p ^ e"'*'  ^ . 


(E-31) 


To  determine  the  variation  of  density  7 under  the  steady  condition 
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note  that  the  use  of  (E-22)  enables  us  to  change  3/9j  to  3/3^;  we  evaluate  ? in  (9)  by  use  of  {E-26) 
and  (E-31);  then,  (12)  yields 

Therefore,  if  the  density  in  the  active  region  oab  (Fig.  6a)  is  constant,  denoted  by  yg,  the  dilatation 
occurs  as  described  by  the  equation 

7 = Yo  e~^  ^ . 

If  the  dilatation  continues  at  the  rate  given  by  this  equation,  the  density  7p  in  the  passive  region 
must  be 

. (E-32) 

Ifp  = 30°,  this  equation  gives  7p/7o  = 0.402.  If  the  dilatation  is  too  large,  the  minimum  density 
7,  must  be  reached  at  angle  i/zq, 

n = (E-33) 

and  maintained  in  the  region  i/z  ^ ^0- 

Point  0 is  a singular  point,  at  which  (E-31 ) yields 

0=0°. 


We  cannot  yet  discuss  the  mathematical  modification  necessitated  by  the  physical  impossibility  of 
the  singularity  at  point  0. 

Velocity  equations  in  the  constant  density  region 

If  the  flow  is  steady  and  e = 0,  density  y is  constant,  as  (1 2)  shows  us,  along  the  flow  lines. 
Therefore,  x in  the  +m-coincident  constant  density  region  must  be  given  by 

X = p/2  (E-34) 

where  (9)  is  used.  Therefore,  the  e<haracteristic  curves  in  this  case  are  orthogonal. 

The  +m-coincident  e-singlets  are  given  by 

^ = cot(i//-p-6) 

where  (1 8),  (E-1 ),  and  (E-34)  are  used.  Eliminating  i/z  by  use  of  (E-2),  this  equation  becomes  a 
total  differential  equation 

{ydx-x  dy)  + (x  </x  +y  dy)  tan  p = 0. 

This  integrates  to 

r = r\  e'** 

by  use  of  (E-4),  where  tj  is  the  initial  value  of  r in  the  constant-density  region.  In  this  case  we  express 
the  doublet  in  (E-5)  with 
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where  ^ is  the  value  of  r at  = i//o-  Solving  the  above  two  equations  we  find 

= X log(T}/|)  (E-35) 


where 


X = sin/ocosp. 


VVe  introduce  the  curvilinear  coordinates  | and  tj  in  place  of  the  lengths  s and  t of  the  doublet 
and  the  e-singlet,  respectively.  We  transform  (E-12)  and  (E-14)  to 

■» 

where  {E-34)  is  used.  These  equations  reduce  to 

dv, 


dy 

dv2 


+ v,  = 0 


+ v,  = 0 


dx 

when  we  put 

X = X log^ 
y = XIogT?. 


(E-36) 

(E-37) 


(E-38) 


Thus,  we  find  that  and  are  governed  by  the  equations 
9^  V, 

— -J.  — I/. 

9x9y  * 

9^  v'j 

" "'2- 


(E-39) 

(E-40) 


These  are  the  same  types  of  hyperbolic  differential  equations  whose  boundary-value  problems  are 
completely  solved  (Courant  and  Hilbert,^  p.  458). 


Solution  in  the  first  constant-density  subregion 

Let  the  line  of  in  Figure  7 be  the  initial  line  at  which  the  constant  density  region  begins.  Draw 
the  e-singlet  passing  through  f ; and  let  g be  the  intersection  of  the  e'-singlet  with  the  line  oc.  The 
velocity  and  its  derivatives  must  be  continuous  across  of,  because  this  line  is  noncharacteristic  in 
the  velocity  field.  The  unique  solution  of  (E-39)  must  be  found  by  use  of  the  boundary  conditions 


V 


1 


I'. 


9i') 

9x 


= 0, 


and 


9v 

ly 
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2nd 

Subregion 


Figure  7.  The  first  and  second  subregions  of  the  constant  density  region. 


on  of,  where  1(1  = ij/Q,  i.e.,  x = y. 
Assume  the  solution 

V,  = F(t) 


t = y-x. 

Then  (E-39)  becomes 
dt^ 

Thus  we  find  the  solution 


V]  = Vcos{y-x). 

Use  of  (E-36)  yields 

i'2  = V'sin(y-x). 

Use  of  (E-38)  and  (E-35)  transforms  them  to 
= Vco%(\l/ -iI/q). 
i'2  = Vs\n(\l)~\l)Q). 

Thus,  die  velocity  vector  v is  given  by 

V = V[-  c^  sin(^o  S)  - 'y  cos(^o  ®)1 

where  (E-10),  (E-7),  and  {E-8)  are  used.  The  direction  of  v,  therefore,  is  u,(,  evaluated  at  i//  = \I)q. 
Therefore,  the  velocity  in  the  region  ofg  is  constant,  being  equal  to  the  initial  velocity  on  the  line  of. 

Solution  in  the  second  constant  density  subregion 

The  solution  in  the  subregion  ^ needs  the  use  of  the  Riemann  function 

= Jo  WiF -P)(<l-y)  1 (E-41 

which  satisfies 


J 


f 


= w 

dxdy 

where  /q(  ) is  the  Bessel  function  of  zeroth  order. 

Use  of  the  x,  y coordinates  defined  by  (E-38)  maps 
the  subregion  ^ in  Figure  7 to  the  region  FGC  in 
Figure  8.  In  more  detail,  fc(|  = C)  is  mapped  to 
FC(x  = a);fi(r}  = 8)  to  FG(y  = (7);andgct^  = (7r/2)] 
to  GCLv  - X = (ff/2)  - ^ol  > "'here  fi  is  the  length  of 
of  and  (7  = X log  8.  In  the  following,  we  shall  solve 
(E-39)  on  the  assumption  that  v^and  its  normal  deriva- 
tive across  the  characteristic  boundary  fg  are  continuous. 
This  assumption  holds  true  for  the  normal  component 
V,  but  does  not  for  the  tangential  component  Vj  on  the 
characteristic  boundary  fg. 

Riemann’s  method  in  this  case  consists  in  the  trans- 
formation of  the  identity 


Figure  8.  Solution  in  the  second  sub- 
region. 


If 

□ PQRF 


dxdy  ^ dxdy 


dxdy  =0 


where  the  range  of  integration  is  the  rectangle  PQRF  in 
Figure  8.  Point  P is  an  arbitrary  point  with  coordinates 
(p,  q),  and  PQ  and  PR  are  parallel  to y-  and  x-axes,  respectively.  The  above  equation  transforms  to 


pnoc  \ ’ ^y) 


dxdy  = 0 


□PQRF I 
which,  on  integration,  becomes 


9tv\ 

3>'/x=p 


dy  = 0. 


Because 


^(x,q)  = 1,  v^(a,y)  = I/,  and  ^ (p,y)  = 0 


we  can  integrate  once  more  to  find 


V,  (p,  q)  - Vw(a,  o)-  J w{x,  a)  (x,  a)  dx 


where  |3w/3y)(p,  y)  and  (3i/,  /dx)(x,  a)  are  the  notations  meaning  to  let  x = p and y = o in 
(9/3y)w(x,  y)  and  (3/3x)v,  (x,  y),  respectively.  Changing  the  dummy  variables,  we  find  the  solution 
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•ffissapi 


ix,y)  = Jo['^'J(<i-x){y-a)]  - J sin(<7  - 1)  Jo[2\J(t  - x){y  - a)  ] dt. 

X 

Use  of  (E-36)  yields 

■y  ‘'2  y)  = -/l  (2v/(o-x)ty-o)  ] - 

a 

J \Jt-x  S)n(ff  - /)  y , [2y/(f  - J()iy  -o)]dt. 


(E-42) 


(E-43) 


The  boundary  condition  ^2  = 0 on  fc  is  satisfied.  One  can  prove  that  V2  is  continuous  on  fg  but 
the  normal  derivative  dvj/Bx  js  not. 

Solution  in  the  passive  region 

The  density  on  a moving  radial  decreases  as  it  moves  from  ob  to  oc  in  Figure  9.  If  7p  in  (E-32) 
is  larger  or  equal  to  Ti,  Tp  is  actually  reached  and  the  constant-density  region  disappears.  In  this 
case,  the  velocity  in  the  passive  region  ocd  is  the  simplest  solution  (E-25),  where  v,  and  ^2  the 
velocity  components  along  the  doublet  direction  and  the  normal  to  the  doublet  direction,  respectively. 

If  the  minimum  density  7^  is  larger  than  7p  in  (E-32),  7,  is  reached,  say  at  of  in  Figure  9,  and 
7p  cannot  be  reached,  then,  the  velocity  in  the  passive  region  becomes  complicated,  as  explained  in 
the  following. 

In  the  passive  region,  the  doublets  and  the  e-singlets  are  straight  lines;  and  velocity  components 
along  these  lines  are  constant.  Therefore,  © = 0 in  this  case;  and  we  have  no  means  of  determining 
X.  However,  for  the  sake  of  presentation,  we  just  assume  x = p/2  and  describe  the  following  velocity 
field. 

The  velocj^  field  in  the  passive  region  is  determined  by  the  velocity  components  on  the  bound- 
aries oc  and  cd  (Fig.  9).  It  is  divided  into  various  segments.  We  draw  the  doublet  straight  lines  and 
the  eaglet  straight  lines  (orthogonal  to  the  former)  at  the  point  g,  and  denote  their  intersections 
with  od  by  points  j_and  h,  respectively.  We  draw  ck  and  jm  parallel  to  the  e-singlet,  and  denote  their 
intersections  with  od  and  oc  by  points  k and  m,  respectively.  Finally,  we  draw  parallel  to  the 
doublet  and  denote  the  intersection  with  oc  by  point  n.  The  velocity  components  on  oc,  cd,  and 
od  are  shown  in  Table  I. 
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Table  I.  Velocity  components  on  the  boundaries  of  the  passive  refion. 

In  this  table  vj  (gc),  i'2  (gt).  and  tbe  similar  symbols  with  other  segment  names 
enclosed  in  the  brackets  are  the  solutions  (£-42)  and  (£-43)  evaluated  on  the 
respective  segments. 


‘'1 

^2 

E3 

V 

0 

og 

ysin  i)/Q 

V cos*o 

gf 

Vi  (gc) 

^2  (gc) 

SK 

V sin  <iiQ 

l^cosij/o 

BT 

V sin  tliQ 

V2(gm) 

ilT 

vj  (gn) 

vjlmc) 

V,  (nc) 

0 

CONCLUSION 

The  exposition  in  this  paper  reveals  that  the  formulation  of  the  in-plane  deformation  of  the  non- 
coaxial Coulomb  material  has  reached  the  stage  of  completion.  However,  it  is  disclosed  that  the 
assignment  of  a set  of  boundary  conditions  cannot  lead  us  to  the  unique  determination  of  the  non- 
coaxial angle;  in  addition,  the  mode  of  flow  must  be  given  to  determine  the  non<oaxial  angle 
uniquely. 

We  can  propose  several  courses  of  further  improvement,  as  explained  below,  that  may  remedy 
this  defect.  However,  the  non-coaxial  angle  is  contained  in  the  constitutive  equations  at  any  stage 
of  the  proposed  improvement;  therefore,  the  constitutive  equations  of  the  Coulomb  material  are 
more  than  material  constant,  even  in  the  final  satisfactory,  if  possible,  formulation. 

Two  courses  suggested  by  the  rational  mechanics  formulation  are  looming  in  our  prospective  re- 
search area.  The  first  is  the  introduction  of  couple  stress.  Mandl  and  Luque'*  interpret  the  non- 
coaxial angle  to  be  induced  by  rotational  movement  of  soil  particles.  It  is  true  that  the  deformation 
of  any  granular  material  is  characterized  by  independent  rotation  of  each  particle.  Then,  it  is 
imperative  for  us  to  introduce  couple  stress  into  the  mechanics  of  granular  material.  So  far,  how- 
ever, couple  stress  has  been  introduced  only  in  ealstic  material  (Mindlin,**  Toupin**).  It  is  known 
that  the  tangential  stress  in  a continuum  endowed  with  couple  stress  is  not  symmetric.  We  believe 
that  the  mathematical  methods  — the  theory  of  characteristic  lines,  the  dyadic  notation  of  tensors, 
and  the  geometric  expression  of  tensors  — that  have  been  successfully  used  in  this  paper  can  be  ex- 
tended to  the  nonsymmetric  system. 

The  second  course  is  the  introduction  of  the  pore  space  distribution.  Goodman  and  Cowin’ 
developed  a mathematical  concept  of  pore  space  distribution  and  introduced  it  into  the  frame  of 
elasticity.  In  the  contemporary  soil  mechanics,  no  rational  method  is  employed  to  conceptually 
convert  a collection  of  particles  to  a continuum.  Goodman  and  Cowin’s  approach  supplies  an  ideal 
correction  to  the  traditional  defect.  Snow  mechanics,  which  is  characterized  by  extraordinarily  large 
pore  space,  especially  needs  this  correction.  However,  we  do  not  know  yet  how  plastic  constitutive 
equations  can  be  introduced  in  this  new  system. 

Mandl  and  Luque,'*  by  citing  Truesdell,“  showed  that  in  a general  isotropic  continuum  the  non- 
coaxiality of  the  principal  axes  of  stress  tensor  and  strain  rate  tensor  occurs  only  during  in-plane, 
but  never  during  three-dimensional,  deformation.  In  other  words,  to  discuss  the  three-dimensional 
plastic  deformation  of  the  non<oaxial  Coulomb  material,  we  must  determine  a plane  in  which  the 
in-plane  deformation  occurs.  The  in-plane  deformation  plane  may  develop  fully  or  partially.  The 
complexity  of  the  deformation  of  the  model  we  are  proposing  is,  we  hope,  the  exact  reflection  of 
the  complexity  of  the  actual  three-dimensional  deformation. 
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The  three-dimensional  yield  criterion  of  the  Coulomb  material  was  given  in  Takagi^'  by  applying 
the  theory  of  plastic  potential,  as  a one-parameter  continuum  connecting  Tresca's  yield  aiterion  and 
von  Mises’  yield  criterion.  It  may  be  possible  that,  rather  than  giving  a specific  value  to  the  parameter 
and  picking  up  a specific  three-dimensional  yield  criterion,  the  parameter  of  this  continuum  may  be 
an  additional  unknown  for  describing  the  three-dimensional  deformation  of  the  Coulomb  material. 
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